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SECTION  1 


INTRODUCTION 


On  30  October  1962  the  United  States  conducted  a  nuclear  test  of  a  near¬ 
surface,  high-yield  device,  called  Housatonic,  which  was  detonated  at  1600  GMT  above 
the  Pacific  Ocean  near  Johnston  Island.  Observations  made  during  this  test  suggested 
that  the  explosion  generated  a  traveling  ionospheric  disturbance,  which  has  been  inter¬ 
preted  as  an  acoustic-gravity  wave  (AGW)  (see  reference  4).  This  interpretation  has 
resulted,  for  the  most  part,  from  the  analysis  of  critical  frequency  and  virtual  height 
measurements  of  vertical-incidence  ionograms  made  at  various  ionosonde  stations  in 
the  surrounding  Pacific  region.  Similar  observations  made  following  the  58  MT  Soviet 
test  at  Novaya  Zemyla  on  30  October  1961  also  detected  ionospheric  disturbances  over 
Europe  later  interpreted  as  produced  by  an  acoustic-gravity  wave  emanating  from  the 
explosion. 


Acoustic-gravity  waves  in  the  upper  atmosphere  can  produce  variations  in 
ionospheric  electron  densities  great  enough  to  cause  changes  in  the  path  of  a  high 
frequency  (HF)  signal  propagating  through  the  ionosphere.  These  changes,  if  severe 
enough,  may  disrupt  HF  communications  in  the  affected  regions.  For  the  two  cases 
mentioned  above,  critical  frequency  variations  were  as  large  as  50%  and  virtual  height 
variations  as  large  as  100  km.  Such  variations  are  not  too  severe  in  general  to  cause 
a  blackout  in  HF  communications.  However,  if  many  more  similar  bursts  occur  in  a 
given  region,  and  sufficiently  closely  spaced  in  time,  then  ionospheric  electron  density 
changes  may  be  severe  enough  to  cause  a  major  disruption  in  HF  communications. 
Furthermore,  if  enough  acoustic-gravity  waves  produced  by  many  low-altitude  nuclear 
explosions  interact  in  some  region,  spatial  variations  in  ionospheric  electron  density 
with  length  scales  on  the  order  of  10  km  may  develop,  which  may  produce  enough 
scatter  in  an  HF  radio  wave  to  cause  am  unacceptable  degradation  in  signal  strength. 

This  report  describes  a  model  which  computes  ionospheric  electron  density 
changes  produced  by  acoustic-gravity  waves  generated  by  any  number  of  high-yield, 
low-altitude  nuclear  explosions.  The  main  requirement,  of  this  model  is  that  it  con¬ 
tain  only  fast-running  analytic  formulae  possessing  second  order  derivatives  in  space 
and  time,  which  can  also  be  determined  and  coded  as  analytic  formulas.  The  im¬ 
plementation  of  this  model  in  the  RAYTRACE  code,  which  must  not  be  appreciably 
slowed-down  by  the  AGW  model,  ma&es  this  requirement  an  absolute  necessity  (see 
the  companion  report,  reference  6,  for  details  on  the  implementation  and  use  of  the 
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AGW  model  in  ray  tracing  studies).  The  model  is  based  on  the  linearized  hydrody¬ 
namic  theory  of  acoustic-gravity  waves  (see  reference  2).  This  theory  is  used  to  justify 
a  single  burst  model  that  fits  well  with  the  data  obtained  for  the  two  events  men¬ 
tioned  above.  Section  2  provides  a  short  description  of  the  data,  concentrating  on 
those  aspects  of  the  data  which  are  most  important  to  the  modeling  process,  while 
section  3  presents  a  thorough  outline  of  the  theory  of  acoustic-gravity  waves.  Section  4 
then  shows  how  this  theory  leads  to  predictions  about  hydrodynamic  density  changes 
produced  in  the  upper  atmosphere  by  an  AGW  emanating  from  a  localized  source.  A 
comparison  with  the  data  is  made  at  this  point. 

Using  an  entirely  linear  approach  to  compute  density  changes  eventually 
leads  to  difficulty  if  either  the  size  of  an  individual  burst  is  too  large  or  if  the  effects 
from  a  large  number  of  interacting  gravity  waves  are  attempted  to  be  computed  by 
linear  superposition  of  the  single  burst  case.  For  a  multiburst  model,  an  approach 
must  be  taken  that  uses  as  much  as  possible  of  the  single  burst  model  but  which 
gives  reasonable  results  (e.g.  positive  densities)  when  the  effects  are  combined.  These 
issues  are  resolved  in  section  5,  where  a  method  of  combining  single  burst  effects  is 
described. 

As  the  number  of  bursts  increases  and  ionospheric  electron  motions  be¬ 
come  correspondingly  large,  variations  in  electron  density  occur  not  only  from  volume 
changes  produced  by  hydrodynamic  motions,  but  also  occur  by  diffusion,  chemical  in¬ 
teraction  with  the  atmosphere,  as  well  as  solar  production  during  daytime.  Section  6 
explains  how  these  ideas  are  incorporated  into  the  model.  Finally,  section  7  presents 
the  final  version  of  the  model,  employing  the  ideas  described  in  the  earlier  sections. 
Some  examples  are  provided  to  demonstrate  the  type  of  results  to  be  expected  from 
using  the  model. 
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SECTION  2 


REVIEW  OF  DATA 


The  Housatonic  device  was  detonated  near  Johnston  Island  at  approximately 
1600  GMT  (4:00  a.m.  local  time).  Variations  in  critical  frequency  of  the  F2  layer 
of  the  ionosphere,  foF2,  were  observed  at  a  number  of  sounding  stations  employing 
vertical-incidence  ionosondes  at  various  locations  in  the  Pacific  Ocean.  Virtual  height 
variations  were  also  observed.  These  variations  had  the  form  of  oscillations  in  the 
measured  quantity  about  its  ambient  value,  with  a  period  of  oscillation  suggesting  the 
passing  of  an  acoustic-gravity  wave.  The  difference  in  the  time  between  the  onset  of 
these  oscillations  and  the  time  of  the  burst,  which  is  related  to  the  speed  at  which 
a  disturbance  generated  by  the  explosion  travels  to  the  observation  point,  are  also 
consistent  with  a  gravity  wave  interpretation.  Details  supporting  these  remarks  are 
given  below. 

The  earth’s  ionosphere  contains  free  electrons  with  densities  that  vary  with 
altitude.  The  ionospheric  electron  density  generally  increases  with  altitude  until  a 
maximum  is  reached  in  the  F  region,  located  about  300  km  above  the  earth’s  surface. 
Vertical-incidence  ionosondes  can  determine  this  maximum  density,  also  called  the 
peak  density,  by  sending  an  electromagnetic  signal,  whose  frequency  is  continuously 
swept  from  about  1  to  20  MHz  (in  the  HF  range),  vertically  upward  into  the  ionosphere 
and  measuring  the  time  it  takes  the  signal  to  return.  The  maximum  frequency  for 
which  there  is  a  return  signal  is  called  the  critical  frequency,  foF2.  (During  daytime 
the  F  region  usually  has  two  regions  with  local  maxima,  called  the  FI  and  F2  layers. 
The  upper  layer,  which  has  the  greatest  density  and  remains  during  nighttime,  is  the 
F2  layer.  The  critical  frequency,  foF2,  measures  this  maximum  density.)  The  critical, 
or  plasma,  frequency  of  the  F2  layer,  foF2,  is  related  to  the  peak  electron  density, 
ne,m«LX>  hy  the  equation 


foF2[in  MHz]  =  2.84  \/nCimax(in  units  of  10s  cm-3].  (1) 

The  virtual  height,  h1,  of  electrons  of  a  given  plasma  frequency,  is  given  by 
h1  =  cf/2,  where  c  is  the  speed  of  light  in  free  space  and  t  is  the  time  it  takes  the  HF 
signal  of  frequency  equal  to  the  plasma  frequency  to  reflect  from  those  electrons  and 
return  to  the  ground.  This  is  not  equal  to  the  true  height  of  the  electrons  because 
the  HF  signal  slows  down  as  it  travels  through  higher  electron  densities.  The  index  of 
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refraction,  n,  of  an  electromagnetic  wave  of  frequency,  /,  in  a  plasma  (in  the  absence 
of  collisions  and  magnetic  field)  with  plasma  frequency,  fp,  is 


n  = 


(2) 


It  can  be  seen  that  as  a  wave  of  given  frequency  (/  >  /p)  encounters  plasma  of  greater 
density,  and  hence  greater  plasma  frequency,  the  index  of  refraction  approaches  zero, 
and  the  wave  will  reflect  when  /  =  fp.  The  true  height,  h,  is  related  to  the  virtual 
height  by  the  integral 


ti  = 


dz 

\J  1  -  ne(z)/ne{h) 


(3) 


where  ne(z)  is  the  electron  density  as  a  function  of  altitude,  z 

Figure  1  shows  a  plot  of  the  critical  frequency,  foF2,  as  a  function  of  time 
(GMT)  above  the  sounding  station  located  on  the  island  of  Tern,  about  1200  km  north 
of  ground  zero  for  Housatonic.  The  solid  line  represents  the  data  while  the  dotted  line 
shows  the  usual,  ambient  value  of  foF2  for  those  times  of  day.  The  ambient  value  is 
rising  sharply  at  these  times  due  to  sunrise,  which  occurred  at  about  1750  GMT  (5:30 
a.m.  local  time).  It  is  readily  observed  that  the  measured  values  of  foF2  oscillate 
about  the  ambient  value.  The  time  of  the  detonation,  labelled  T0,  is  1600  GMT, 
while  the  first  occurrence  of  the  oscillation  is  about  30  minutes  later.  This  travel  time 
corresponds  to  an  average  speed  of  about  .6  km/sec,  which  is  near  the  average  speed 
of  sound  along  a  line  from  ground  zero  to  an  observation  point  in  the  F  region  above 
Tern.  Notice  that  the  number  of  oscillations  is  about  two  or  three  before  damping 
out  and  that  the  time  elapsed  between  the  first  crest  and  the  second  crest,  labelled  T 
in  the  figure,  is  about  50  minutes.  Also  notice  that  the  second  oscillation  appears  to 
be  broader  than  the  first,  suggesting  that  the  period  of  oscillation  is  increasing  with 
time. 


Although  not  shown  in  the  figure,  the  virtual  height  also  showed  some  oscil¬ 
lation,  with  a  decrease  in  virtual  height  corresponding  to  the  first  increase  in  critical 
frequency.  As  will  be  seen,  the  orientation  of  the  earth’s  magnetic  field  lines  relative  to 
the  direction  of  propagation  of  the  gravity  wave  determines  whether  electrons  travel 
up  or  down  the  field  lines  as  the  gravity  wave  passes.  For  Tern,  the  field  lines  are 
oriented  such  that  electron  motion  down  the  field  is  expected  to  occur  as  the  gravity 
wave  first  passes,  as  is  observed  (see  reference  4). 
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Figure  2  shows  a  similar  plot  of  the  critical  frequency,  foF2,  as  a  function  of 
time  (GMT)  above  the  sounding  station  located  on  the  island  of  Tonga,  about  3900 
km  north  of  ground  zero.  This  plot  also  shows  an  oscillation  of  foF2  about  the  ambient 
value.  The  onset  time  of  the  oscillation  is  about  90  minutes  following  the  burst,  which 
gives  an  average  travel  speed  of  about  .7  km/sec,  slightly  larger  than  the  value  for 
Tern.  The  number  of  oscillations  observed  is  also  two  or  three.  The  time  between 
the  first  crest  and  the  second  crest,  labelled  T,  is  about  100  minutes,  greater  than 
the  corresponding  value  for  Tern.  As  in  the  Tern  case,  the  time  between  oscillations 
appears  to  increase  with  time.  Finally,  the  relative  change  in  foF2  with  respect  to 
ambient,  as  measured  by  the  amplitude  of  the  first  oscillation,  is  smaller  for  the  Tonga 
data  than  for  the  Tern  data.  As  for  Tern,  the  virtual  height  decreases  initially,  as  is 
expected  from  the  orientation  of  the  earth’s  magnetic  field  above  Tonga.  Although 
these  two  sites  provided  the  cleanest  data,  similar  variations  in  foF2  were  observed  at 
other  locations,  the  details  of  which  can  be  obtained  in  reference  4. 

Figure  3  shows  plots  of  the  critical  frequency,  foF2,  as  a  function  of  time 
(MET)  above  various  European  sounding  stations,  following  the  58  MT  Soviet  test 
at  Novaya  Zemyla  on  30  October  1961.  The  stations  shown  in  the  figure  are  located 
at  increasing  distance  from  ground  zero  going  up  in  the  figure,  from  about  1400 
km  for  Kiruna  to  about  4300  km  for  Athens.  Although  the  data  used  to  make  the 
figure  were  not  very  good,  some  general  observations  can  be  made.  The  most  obvious 
characteristic  is  that  the  period  of  oscillation  increases  with  increasing  distance  from 
ground  zero.  It  is  also  seen  that  the  number  of  oscillations  at  most  locations  is  about 
four  or  five.  A  careful  examination  also  shows  that  there  was  an  initial  decrease  in 
foF2,  with  a  corresponding  increase  in  virtual  height  (not  shown  in  the  figure),  which 
is  expected  from  the  magnetic  field  orientation  at  those  sites.  It  is  not  clear  from  the 
figure  whether  or  not  there  is  much  decrease  in  the  amplitude  of  the  oscillations  with 
increasing  distance  from  ground  zero. 

The  data  shows  a  number  of  features  which  will  be  important  in  constructing 
an  acoustic-gravity  wave  model.  The  most  important  general  feature  is  that  the 
F  region  of  the  ionosphere  oscillates  at  quite  distant  points  following  a  large,  near- 
surface  explosion.  The  period  of  oscillation  increases  with  increasing  distance  from  the 
source,  approximately  linearly  with  distance.  The  period  of  oscillation  at  any  given 
point  also  increases  slightly  with  time.  The  number  of  oscillations  is  approximately 
constant  for  a  given  burst,  independent  of  location.  The  amplitude  of  the  oscillat'on 
decreases  with  distance,  although  the  data  are  insufficient  to  determine  what  the 
scaling  law  with  distance  is.  The  direction  of  ionospheric  electron  motion  depends 
on  the  orientation  of  the  earth’s  magnetic  field  lines  with  respect  to  the  direction  of 
propagation  of  the  gravity  wave.  In  section  4  it  will  be  shown  how  the  general  theory 
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of  acoustic-gravity  waves  can  be  used  to  derive  these  properties,  which  will  be  the 
basis  for  constructing  the  AGW  model.  Before  presenting  these  derivations,  a  review 
of  the  general  theory  of  acoustic-gravity  waves  will  be  presented  in  section  3. 
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SECTION  3 


THEORY  OF  ACOUSTIC-GRAVITY  WAVES 


Most  properties  of  acoustic-gravity  waves,  in  particular  those  discussed  in 
the  previous  chapter,  can  be  derived  by  assuming  that  the  earth’s  atmosphere  can 
be  described  as  a  stationary,  stratified  fluid,  with  constant  temperature  in  a  uniform 
gravitational  field.  Acoustic-gravity  waves  come  about  as  small  (first  order)  pertur¬ 
bations  about  isothermal  equilibrium  of  the  equations  of  hydrodynamics.  The  fluid  is 
described  by  the  pressure,  p,  density  p,  and  velocity  v,  given  as  functions  of  position 
and  time.  The  equations  of  motion  are  the  equation  of  continuity,  balance  of  linear 
momentum,  and  adiabatic  equation,  given,  respectively,  by 


Dp 

~  +  pV  -v  =  0, 

(4) 

D\ 

pm+*p  =  ps • 

(5) 

Dp 

— +  T,pVv  =  0, 

(6) 

where  the  convective  derivative  is  given  by  D/Dt  =  d/dt  +  v  •  V,  -y  is  the  adiabatic 
constant,  and  gravity  is  assumed  to  be  constant  in  the  negative  z  direction,  g  =  —gex. 
Equation  (6)  is  derived  by  assuming  the  specific  entropy  of  a  material  particle,  s,  is 
constant, 


Ds 

Dt 


=  0, 


(?) 


and  that  the  fluid  is  ideal  (or  perfect), 

P 

s  =  c„ln— , 
P1 


(8) 


and  using  the  continuity  equation  (4). 
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The  density  of  a  fluid  in  static,  motionless  (v  =0),  isothermal  equilibrium 
is  distributed  exponentially  as  a  function  of  height, 


Poe 


-m/H 


where  p0  is  the  density  at  z  =  0.  The  corresponding  pressure  is 


(9) 


P(0>  = 


(10) 


According  to  equation  (5)  the  constants  p0  and  po  must  be  related  to  the  scale  height, 
H,  by 


Po  =  p0gH. 


(11) 


The  speed  of  sound,  c,  is  given  by 

c2  =  T-  =  7—  =  'igH.  (12) 

P  Po 

To  derive  the  linearized  equations  of  motion  from  equations  (4)-(6),  which 
describe  small  perturbations  about  the  solutions  given  by  equations  (11)  and  (12) 
(with  v(°l  =  0),  we  will  use  perturbation  theory  in  a  systematic  way  so  the  interested 
reader  may  use  this  approach  as  a  starting  point  for  higher  order  corrections,  which 
will  not  be  pursued  here.  It  is  easiest  to  work  in  cartesian  coordinates  and  suppress 
the  y  component  of  the  velocity,  which  can  be  included  at  any  point  in  what  follows 
by  making  the  obvious  changes.  Therefore,  let  the  velocity,  v,  be  written 


v(x,z,t)  =  u(x,z,t)ex  +  w(x,z,t)ex. 
The  equations  of  motion  in  component  form  are 


(13) 


dp  dp  du  dw 

dx  dz^  Pdx  P  dz  ’ 

(14) 

du  du  1  dp 

+  — h  W— — b  =  0, 

dx  dz  o  dx 

(15) 
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(16) 


dw  dw 

~dt+UJi 


dw  1  dp 

4-  W -  -j - 

dz  Pdz 


■ 9 , 


dp 

dt 


t— — r  w— - 
3x  dz 


(17) 


A  solution  to  these  equations  is  assumed  in  the  form  of  an  expansion  in 
powers  of  the  perturbation  parameter,  s, 


u(x,z,t)  =  eu  (1)(x,  z,  t)  +  z  VJ)(x,  z, t)  H - , 


(18) 


u/(x,  z,  t)  =  zto  ^(x,  z,  f )  +  zlu;^(x, z,t)  H - ,  (19) 

p(x,  z,  f)  =  p(0)  4-  ep{1)(x,  z,  t)  +  e2p{-2\x,  z,  t)  +  •  • « ,  (20) 


p(x,  z,  t)  =  p(0)  +  ep(1)(x,z,t)  +  zVJ)(x,z,t)  +  •••,  (21) 

where  p^  and  p(°)  are  given  in  equations  (9)  and  (10).  Also  needed  is  the  expansion 
for  the  function  p~l, 


Substitution  of  equations  (18)  to  (22)  into  equation  (14)  and  using  the  con¬ 
densed  notation,  uz  =  du/dx ,  gives 


ePt(1)  +  e2p|2)  H - +  +  •  •  ■)  (epW  -| - )  + 

(eu/W  +  e2wW  -\ - )  +  ep^  +  •••)  + 

(p(°)  +  epW  -| - j  (eu^  +  ew +  e2u^  +  zJu/^  -f - j  =0.  (23) 

Collecting  terms  in  like  powers  of  e  gives 
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e  (pj1^  +  v/WpW  +  p{0>ttpJ  4-  p(°)«;^)  4 

eJ  (pjZ^  4-  uWpW«Wpt°J  +  w>WpW  4-  p(°>t4J)  4-  p<°>ti>«  4-  p^u^  4-  p^tup^ 

4-  ■  •  •  =  0.  (24) 

Substitution  of  equations  (18)  to  (22)  into  equations  (15)-(17)  produces  three  more 
equations  similar  to  equation  (24).  Setting  the  coefficient  of  each  power  of  e  equal  to 
zero  and  assembling  the  sets  of  equations  that  result  order  by  order  yields  through 
second  order, 

Zeroth  order: 


-p(0> 


9 , 


from  which  it  can  be  concluded  that  po  =  PoQH  as  in  equation  (11). 
First  order: 


(25) 


pj1*  4-  4-  P(0)up>  4-  p*°>»W  =  0, 

(26) 

4>'  +  ^p<'>  =  o. 

(27) 

W*  )  +  p(0)p(0)Pi0)  =  °> 

(28) 

Pt^  4-  w^pp*  4-  7P^up^  4-  'ip^w^  =  0, 

(29) 

which  are  the  usual  linear  AGW  equations. 
Second  order: 


p|2^  +  tv^pW  +  p(°)tt£a)  4-  p<°>«4J>  =  -v^pP*  -  ti/^pP*  -  p^up)  4-  p^typ),  (30) 
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(31) 


«P> + =  -»<*>«<■>  -  «-(V>  +  ^ 


W*  +  „(°)P*  «{0)„(0) 


,J_£!!L(°)  =  -uUWiJ-^W^W+JL^pti) _ L  ( JO)  f32l 

p(0)  p(°)Pt  *  '  y<V0)P’  p(0Hp(0)/  P*  ’ 


+  vMpW  +  <7p^«^  +  =  —u^p\P  —  w^pM  —  -yp(1)uW  4-  'ip^wW,  (33) 


which  are  also  linear  equations  in  the  second  order  quantities  appearing  on  the  left 
hand  sides.  The  functions  on  the  right  hand  sides  are  of  no  greater  than  first  order 
and  are  to  be  considered  given  solutions  of  the  first  order  equations.  The  second  order 
equations  are  given  here  for  reference,  they  will  not  be  needed  in  what  follows.  As  a 
check  on  the  expansions,  it  will  be  seen  that  the  sum  of  the  superscripts  of  each  term 
is  equal  to  the  order  of  the  expansion. 

We  will  now  focus  our  attention  on  the  linear  AGW  equations  (26)-(29). 
Solutions  of  these  equations,  subject  to  appropriate  initial  and  boundary  conditions, 
describe  the  motion  of  small  disturbances  from  the  equilibrium  configuration  given  by 
equations  (9)  and  (10).  If  the  disturbance  is  too  large,  then  these  equations  do  not 
provide  an  adequate  description  and  higher  rtrder  corrections  or  an  entirely  different 
approach  are  needed.  It  has  been  found  that  t  m  linear  AGW  equations  usually  provide 
an  adequate  description  of  traveling  ionosphere  disturbances  at  points  distant  from 
the  source  of  the  wave.  In  our  case,  we  wish  to  describe  ionospheric  motion  at  points 
far  from  and  at  greater  altitudes  of  a  point  which  is  the  source  of  a  large  localized 
atmospheric  disturbance.  The  linear  equations  will  not  yield  solutions  at  points  near 
the  source.  Eventually,  however,  the  disturbance  will  have  spread  out  enough  that  an 
initial  configuration  can  be  defined  which  is  close  to  the  linear  regime.  The  difficultly, 
of  course,  is  that  this  initial  configuration  is  unknown.  It  will  be  seen,  however,  that 
much  information  is  revealed  by  making  quite  general  assumptions  about  this  initial 
configuration  and  employing  general  plane  wave  solutions  of  the  linear  equations, 
which  will  now  be  described. 

It  can  be  shown  that  plane  solutions  to  equations  (26)-(29)  exist  if  an  expo¬ 
nential  factor  with  scale  height  2 H  is  added  to  the  velocity  functions  and  exponential 
factors  with  scale  height  —2 H  are  added  to  the  density  and  pressure  changes  (see,  for 
example,  reference  8).  That  is,  there  exist  solutions  of  the  form 
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4"(r, <;  =  u,(k)  c‘/2B  t'(k  ■ r  -  u(kW, 

(34) 

4‘>(r,  ()  =  ui,(k)  «‘(k ' r  “  “W*) , 

(35) 

p£>(r.‘)  =  p,(k)  <‘(k ' r  "  “(kM, 

(36) 

f)  =  Pl(k)  «‘(k  • r  ~  "WO . 

(37) 

Here,  the  complex  constants,  ui,  wlf  p\,  and  pi,  depend  only  on  the  wave  number  k 
and  provide  the  relative  amplitude  and  phase  of  the  first  order  perturbations.  Their 
precise  form  can  be  determined  in  a  straightforward  manner,  but  will  not  be  needed 
here  (see  reference  8  for  details).  The  angular  frequency,  w(k),  is  determined  by  the 
dispersion  relation 


=  +  ?  <38) 

where  kh  stands  for  the  horizontal  wave  number  which  is  kz  for  the  planar  case  where 
there  is  no  y  dependence  or  =  kl  +  k*  for  the  general  or  radially  symmetric  case, 
which  will  be  the  case  in  the  subsequent  presentation.  The  acoustic  cutoff  frequency, 
u>a,  and  the  Brunt-Vaisala  frequency,  w*,  are  given  by 

"2 =(!)’.  ^ =(-»-»)?•  <39) 

The  dispersion  relation  (38)  can  also  be  written 

w*  =  |  (c2*2  +  w*)  ±  {c3k3  +  w2)2  -  4 c2kloj{.  (40) 

It  is  easy  to  show  from  this  equation  that  waves  can  propagate  (w2  >  0)  only  if  u  >  u>„, 
which  is  called  the  acoustic  branch,  or  w  <  ut,  which  is  called  the  gravity  wave  branch. 
In  the  subsequent  analysis,  most  of  our  concern  will  be  with  waves  that  are  composed 
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of  solutions  whose  frequencies  are  in  the  gravity  wave  branch  and  we  will  use  the  term 
gravity  wave  to  describe  such  solutions. 

A  general  solution  to  equations  (26)-(29)  can  be  obtained  as  a  superposition 
of  the  plane  wave  solutions  (34)-(37).  For  example,  the  horizontal  velocity  u(r,f)  can 
be  written 


«(r.<)  =  p^jSTJ  //00“k’(M)& 

If  the  initial  condition  for  u  is  (the  initial  condition  for  the  function  u  is  the  horizontal 
velocity  as  a  function  of  position  following  the  onset  of  the  localized  disturbance  at 
some  later  time,  which  will  be  taken  as  t  =  0,  such  that  the  disturbance  is  small 
enough  so  that  the  linear  approximation  is  valid) 


u(r,0)  =  /( r),  (42) 

then  a  straightforward  application  of  the  Fourier  theorem  shows  that 

/(k) = (4^  4ki  /  ■ r  *•  i43> 

The  initial  conditions  can  be  used  to  determine  the  parameter  e,  usually  taken  to  be 
one. 


With  a  few  simple  assumptions  regarding  the  nature  of  the  function  /( k),  a 
great  deal  can  be  learned  by  studying  equation  (41)  and  using  the  dispersion  relation 
(38),  which  is  the  subject  of  the  next  section. 
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SECTION  4 


ASYMPTOTIC  ANALYSIS  OF  GRAVITY  WAVES 


By  applying  standard  asymptotic  analysis  to  integral  solutions  of  the  AGW 
equations  such  as  that  appearing  in  equation  (41),  a  great  deal  can  be  learned  about 
the  nature  of  the  solutions  at  positions  far  from  a  localized  initial  disturbance  for  late 
times.  It  will  be  seen  that  much  of  the  information  can  be  obtained  from  the  dispersion 
relation  independent  of  the  precise  details  of  the  initial  conditions.  In  particular,  the 
method  of  stationary  phase  will  be  used,  which  can  be  stated  as  follows:  let 


/(x)  =  //(  x)exp[i'A<£(x)]dx,  xaj*!,...,*.).  (44) 

Under  suitable  conditions,  which  can  be  found  in  reference  1,  this  integral  can  be 
expanded  in  an  asymptotic  expansion  for  A  — >  oo,  the  first  term  of  which  is 


/(x o)  exp 

t'A^(xo)  +  i^sgndet  (x0))] 

det  | 

(xo)) 

1/2 

(45) 


where  the  sum  is  over  all  those  points,  Xo,  such  that  V<£(xo)  =  0.  Higher  order  terms 
contain  larger  negative  powers  of  A,  which  approach  zero  faster  than  the  term  above 
as  A  — *  oo. 

Applying  this  well-known  result  to  equation  (41)  provides  an  asymptotic 
approximation  for  the  horizontal  velocity,  and  similarly  for  the  other  hydrodynamic 
field  variables,  as  t  — ►  oo, 


u(r,t) 


/(k0)u1(k0)  exp 

» (k0  • 

r  - 

*;(k0)t)  +  i^sgndet 

f3Mk0n 
(  dkidkj  J 

*3/2 

det 

f^MkoA 
i,  dkidkf  ) 

1/2 

(46) 


where  k0  is  a  solution  of 
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(47) 


duj  _  r*  du>  _  z 

dirh~7’  dk,  ~  7 

with  rj[  =  x*  4-  y2  and  k\  =  fc2  -f  kj. 


4.1  FREQUENCY  OF  OSCILLATION. 

To  solve  equations  (47)  to  determine  w0  =  w(k0)  we  use  the  dispersion 
relation  (38)  to  obtain 


du) 

dk x 


*•*  w 

*r 


kfal 

s 


du 

dk. 


kt 


w 


(48) 


Using  these  relations  in  equation  (47)  and  taking  the  ratio  of  the  second  to  the  first 
equation  gives 


kM 


kh  (u/2  *) 


1 

tan0’ 


(49) 


where  cos#  =  z/r,  rs  =  rj  +  z1.  From  this  and  the  dispersion  relation  it  follows  that 


i  2 


w4  sin2  6 


c2  (u2  —  w2)  (a;2  —  u2  cos 2  ffj 


(50) 


j  __  u>2  -  -  "»)  cos2  6 

e2  (a;*  —  vf  cos2  oj 


(51) 


Taking  the  sum  of  the  squares  of  equations  (47)  gives 
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Using  equations  (50)  and  (51)  in  this  equation  finally  yields 


(52) 


_ ~  tj)  (k*  ~  cos2  6) _  =  _r^_ 

w2  J^w2  —  w2)  (w2  —  w2  cos2  —  (w2  —  w2)  (w2  sin2  j  c  t 

The  solution(s)  of  this  equation  give(s)  those  frequencies,  u>0,  to  be  used  in  the  asymp¬ 
totic  solution  (46),  which  we  will  now  analyze. 

Let  the  square  root  of  the  left  hand  side  of  the  above  equation  be  designated 
/(w)  so  the  equation  is  written 


/(w)  =  (54) 

Solution  of  this  equation  can  be  obtained  in  a  graphical  manner  as  shown  in  figure  4. 
The  solid  curve  in  the  figure  is  a  plot  of  1/ /(w)  as  a  function  of  w/w6.  The  horizontal 
dashed  line  has  an  ordinate  equal  to  t  =  ct/r.  Solutions  of  equation  (53)  are  at  those 
points  where  the  two  curves  intersect.  It  can  be  readily  seen  that  for  t  large  enough 
there  are  three  solutions,  which  will  be  called  wo.i,  w0,2,  and  w0,j.  It  is  obvious  from 
equation  (53),  that  as  t  —*  00, 


Wo.l  —*  We,  Wo, 2  Wo,3  — ►  w0, 


(55) 


where  wc  =  ub  cos  6. 

Low  frequency  gravity  wave  oscillations  whose  frequency  approaches  wc  are 
consistent  with  observations.  For  observation  locations  at  equal  altitudes,  the  fre¬ 
quency,  ue,  decreases  with  distance  from  the  localized  source  as  l/r,  as  is  observed. 
It  also  is  clear  from  figure  4  that  as  t  — ►  00,  wo  — < ►  wc  (from  this  point  we  will  be 
concerned  only  with  wo.i  which  will  be  called  wo)  from  above,  that  is,  the  oscillation 
frequency  decreases  with  time  (the  period  increases)  as  is  also  observed.  To  the  next 
order  in  t  it  can  be  shown  that 
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l/f(w) 


Figure  4.  Graphical  solution  for  asymptotic  frequency. 
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(56) 


where  c  —  (w4 /wa)e,  which  for  a  typical  atmosphere,  is  a  speed  just  slightly  less  than 
c. 


4.2  SCALING  WITH  DISTANCE. 


The  asymptotic  solution  for  the  horizontal  velocity  given  be  equation  (46) 
is  an  oscillatory  function  of  time  whose  frequency  approaches  we  as  t  — *  oo.  Using 
equations  (50)  and  (51)  it  follows  that 


k  -r  = 


r 

(wJ  -  w*)  (wJ  -  w’) 

(wJ  -  w?) 

(57) 


which  approaches  zero  as  t  — ►  oo,  showing  there  is  no  spatial  oscillation  at  late  times. 

The  amplitude  of  the  oscillation  (suppressing  the  exp(z/2H)  factor),  which 
will  be  called  u0,  is  then 


ti0(r,i) 


/(ko)tti(ko) 

f3/2 

det  | 

(  dJw(k0)\ 
l  dkidkj  ) 

1/2 

(58) 


where  we  have  assumed  /  is  spherically  symmetric  and  written  its  argument  as  k0.  A 
bit  of  arithmetic  shows  that  asymptotically  as  t  — *  oo 


which  gives 


z 

Wo  ~  ~  w4. 

r 


(59) 


(60) 
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It  can  also  be  shown  that 


*  ~  -it-, 

\aK3k,j  kKokl0' 


from  which  it  follows  that 


Finally,  «i  ~  1  and,  therefore,  the  scaling  relation  for  uo  is, 


(61) 


(62) 


(63) 


Consider  now  some  point,  call  it  t0l  in  the  oscillation  where  m  cycles  have 
occurred  (m  may  be  a  fraction),  then 


w0fo  =  2irm  or  to  =  — — 

Wi  z 


(64) 


from  which  it  follows  that 


u0(r,t0)  ~  r*  /  *  (65) 

This  equation  shows  that  for  points  at  altitude,  z,  the  number  of  oscillations,  m, 
which  occur  is  the  same  independent  of  the  distance  from  the  localized  source,  r,  if 
/  is  a  function  of  finite  width,  which  is  the  case  for  a  finite  sized  initial  disturbance. 
This  equation  also  shows  that  the  amplitude  of  the  oscillations,  for  a  given  number  of 
oscillations,  scales  with  distance  as  1/r2. 


4.3  TIME  OF  ONSET. 


It  is  well  known  that  wave  energy  travels  at  the  group  velocity, 
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(66) 


du> 


For  acoustic-gravity  waves  this  value  is  bounded  by  the  speed  of  sound,  c.  This  means 
that  if  a  disturbance  is  localized,  then  an  undisturbed  point  will  remain  so  until  a  time 
equal  to  r/c  has  elapsed,  where  r  is  the  distance  between  the  undisturbed  point  and 
the  closest  point  of  the  disturbance.  If  disturbances  are  confined  to  the  gravity  wave 
branch  (u>  <  w»),  then  the  travel  is  speed  is  actually  bounded  by  c  =  (w*/wa)c.  To  see 
this,  note  that  for  horizontal  propagation,  which  is  the  fastest  way  to  cover  a  distance 
R,  the  dispersion  relation  shows  the  horizontal  component  of  the  group  velocity,  ug, 
can  be  written 


u 


—  w 


11 1/S 


(67) 


For  w  <  this  equation  shows  that  ut  has  a  maximum  equal  to  c  for  u?  =  0.  There¬ 
fore,  for  an  initial  disturbance  that  is  localized  at  the  origin  at  t  =  0,  the  soonest  a 
gravity  wave  disturbance  can  reach  a  point  located  at  a  distant  r,  assumed  to  be  large 
compared  to  the  extent  initial  disturbance,  is  about  t  =  r/c. 


4.4  ELECTRON  DENSITY  CHANGES. 


Up  to  this  point  we  have  concentrated  our  attention  on  the  horizontal  com¬ 
ponent  of  the  neutral  velocity  of  the  atmosphere  during  the  passage  of  an  acoustic- 
gravity  wave  (it  can  be  shown  that  the  vertical  velocity  component  is  smaller  by  a 
factor  of  1/r  asymptotically).  Measurable  changes  in  the  ionosphere  are  produced  by 
indirect  coupling  of  the  ionospheric  electrons  to  the  neutral  motion  of  the  atmosphere. 
It  can  be  demonstrated  that  at  F-region  altitudes  the  coupling  is  such  that  the  elec¬ 
tron  velocity  is  equal  to  the  ionic  velocity,  v,,  which  is  equal  to  the  component  of  the 
neutral  velocity  along  the  earth’s  magnetic  field,  i.e. 

v<  =  (v  •  B)  B,  (68) 

where  B  is  a  unit  vector  in  the  direction  of  the  earth’s  magnetic  field.  (Full  justification 
of  this  relation  can  be  found  in  reference  3.)  Let  the  electron  density  change,  to  first 
order,  be  denoted  n'  so  that  the  electron  density  is  given  by 
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MM)  =  Mz)  +  »'(M). 


(69) 


where  no  is  the  ambient  electron  density,  assumed  to  be  stratified.  The  first  order 
continuity  equation  for  electrons  is  then 


^  +  V  •  (nov.)  =  0.  (70) 

For  a  realistic  atmosphere,  dissipation  of  gravity  waves  occur  as  it  propagates 
vertically  upwards,  due  to  ion  drag,  viscosity,  and  heat  conduction.  Attenuation 
of  the  wave  begins  to  occur  at  approximately  F-region  heights  in  such  a  way  that 
the  exp(z/2H)  increase  in  the  neutral  velocity  is  almost  exactly  cancelled  at  greater 
heights  (see  references  5  and  9).  Therefore,  let  the  form  of  the  neutral  velocity  at 
F-region  heights  be 


=  vJ(k‘T-ut). 


v  =  v0e 


(71) 


Assuming  the  same  time  dependence  for  n#,  equation  (70)  gives 


»'  =  --(v-B) 

u} 


(k  •  £)«<>  +  i(eM  ■  B) 


dn0 

dz 


(72) 


At  the  electron  density  peak,  where  dn0/dz  —  0,  it  can  be  seen  that  the  electron 
density  change  is  in  phase  with  the  neutral  velocity,  the  density  increasing  for  elec¬ 
trons  moving  downward  and  decreasing  for  electrons  moving  upward.  Applying  the 
asymptotic  scaling  relations,  equations  (58)-(65),  we  have  for  the  scaling  relation  with 
distance  for  the  amplitude  of  »'  at  the  peak 


«-l/(^),  (73) 

where  K  is  a  constant  that  depends  on  the  magnetic  field  orientation,  and  the 
point  in  the  oscillation  of  concern  (usually  the  peak).  Notice  that  the  amplitude  of 
the  oscillations,  for  a  fixed  peak  altitude,  falls  off  as  1/r,  which  is  consistent  with  the 
data  for  ranges  that  are  not  too  long. 
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4.5  COMPARISON  WITH  DATA. 


At  this  point  we  have  derived  enough  characteristics  of  gravity  waves  em¬ 
anating  from  a  localized  source  to  construct  a  simple  model  that  agrees  well  with 
the  data  presented  in  section  2.  To  do  this  we  must  choose  a  form  for  the  unknown 
function,  /,  which  is  related  to  the  initial  distribution  of  the  hydrodynamic  horizontal 
velocity.  As  explained  above,  the  form  of  the  function  determines  the  number  of  os¬ 
cillations  as  well  as  the  magnitude  of  the  amplitude,  expected  to  be  different  for  each 
initial  disturbance.  The  frequency  of  oscillation,  onset  time,  and  scaling  do  not.  We 
will  try  a  simple  function  of  finite  width,  a  guassian  of  the  form 


with 


f{k)  =  U  exp 


a  = 


2?rm 


(74) 


(75) 


where  the  constants  U  and  m  are  expected  to  depend  on  the  total  energy  of  the  initial 
configuration.  We  will  make  a  fit  to  the  data  to  determine  these  constants. 

In  order  to  account  for  free  electrons  created  by  solar  radiation,  as  seen  in  the 
data  in  figures  1  and  2,  we  simply  add  a  constant  source  term  for  electrons,  Q(t)  to  the 
right  hand  side  of  equation  (72)  which  agrees  with  the  ambient  data  when  no  gravity 
waves  are  present.  With  this  in  mind,  a  simple  numerical  solution  of  equation  (72)  the 
horizontal  neutral  velocity  of  section  4.1,  the  onset  time  determined  in  section  4.3,  and 
the  function  /  given  by  equation  (74)  with  the  constants  suitably  adjusted,  produces 
the  critical  frequency  changes  plotted  in  figures  5  and  6.  The  dashed  lines  correspond 
to  data  plotted  in  figure  1  and  2,  while  the  solid  lines  are  the  results  of  the  model.  As 
can  be  seen  by  in  the  figures,  the  model  reproduces  the  essential  features  of  the  data. 


Comparison  with  the  Soviet  test  data  allows  a  scaling  relation  with  yield  for 
the  amplitude  U  and  the  parameter  m  which  gives  the  number  of  oscillations  to  be 
determined.  It  is  found  that  reasonable  fits  result  from 


U  =  UqY,  m  =  mQW 


(76) 


where  the  values  of  the  constants  are  determined  by  a  best  fit  to  the  data. 
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SECTION  5 


NONLINEAR  ISSUES 


Electron  density  changes  as  computed  by  equations  (70)  or  (72)  are  valid 
only  under  the  approximation  that  the  change  n'  is  small  compared  to  the  ambient 
electron  density  n«.  The  electron  density  changes  depend  on  the  horizontal  neutral 
velocity,  given  asymptotically  in  equation  (46).  The  magnitude  of  this  velocity  is 
determined  implicitly  by  the  function  /  which,  in  turn,  depends  on  the  magnitude 
of  the  initial  condition,  as  given  in  equation  (42).  If  this  velocity  is  large  compared 
with  the  speed  of  sound,  then  the  approximation  which  shows  that  the  solution  of 
the  linear  first  order  equations  (26)-(29)  are  close  to  the  true  solutions  of  equations 
(14)-(17)  is  no  longer  valid.  Large  values  of  the  neutral  velocity  correspond  to  large 
density  changes,  which  will  eventually  lead  to  unphysical  negative  densities  if  the  linear 
formulas  are  used  for  initial  disturbances  of  increasing  magnitude.  (It  is  assumed  that 
the  magnitude  of  the  initial  disturbance  is  some  increasing  function  of  the  energy 
which  generates  the  initial  disturbance.) 

For  a  number  of  localized  disturbances,  which  would  individually  generate 
velocity  distributions  such  as  equation  (46),  the  velocity,  in  the  linear  approximation, 
would  be  given  by  the  sum  of  the  individual  velocities,  and  the  corresponding  electron 
density  would  be  calculated  using  this  sum.  If  this  number  of  initial  disturbances 
becomes  too  large,  however,  then  this  sum  will  eventually  approach,  in  magnitude, 
the  speed  of  sound,  leading  to  similar  difficulties  described  in  the  previous  paragraph. 

To  model  a  large  number  of  localized  disturbances  of  arbitrary  magnitude, 
it  is  necessary  to  deal  with  these  nonlinear  issues.  We  will  not  be  so  ambitious  as  to 
try  to  approximate  corrections  to  the  actual  equations  of  motion,  but  instead  will  tr. 
to  construct  solutions  with  certain  reasonable  properties  which  will  be  stated  below. 
The  results  of  this  section  must,  therefore,  be  looked  upon  in  their  proper  light  as,  at 
best,  educated  guesses.  The  assumptions  which  will  guide  our  approach  to  modeling 
large  disturbances  are: 

1.  The  electron  density  must  always  be  positive. 

2.  The  solutions  must  approach  the  linear  solutions  for  small  initial  disturbances. 

3.  The  magnitude  of  the  velocity  must  approach  a  constant  as  the  altitude  in¬ 
creases. 
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This  last  assumption  is  based  on  remarks  stated  in  section  4.4  that  in  a  realistic 
atmosphere  gravity  waves  become  damped  with  increasing  altitude  in  such  a  way  that 
the  factor  exp(z/2H)  is  almost  exactly  cancelled. 

We  will  now  outline  the  reasoning  which  will  lead  to  a  solution  with  the 
properties  enumerated  above.  The  vertical  velocity  of  electrons,  we,  according  to 
equation  (68)  is 


W«  =  (v.B)(e,.B),  (77) 

where  the  neutral  velocity  is  v  =  ue&  (the  vertical  velocity  component,  weM,  can  be 
added;  it  falls  with  distance  faster  than  tt  by  a  factor  of  1/r).  Remarks  following 
equation  (57)  show  that  asymptotically  the  spatial  phase  of  this  velocity  is  zero  and, 
therefore,  the  spatial  variation  is  mild  except  for  the  exp(z/2H)  factor.  According 
to  third  assumption  above,  this  factor  must  approach  one  for  altitudes  above  the  F- 
region  peak,  which  will  be  called  ze.  We  will  therefore  try  to  find  some  function,  call 
it  g{z),  with  the  property  that  g  — »  1  as  z  — ►  oo  and  g  ~  exp(z/2H)  for  z  <C  ze  and 
write 


t ot  =  g[z)wo. 


(78) 


where 


w0  —  (e*  •  B)(e,  ■  B)e  z/2^u.  (79) 

The  function  w0  is  a  very  mild  function  of  position  and  can  be  considered  a  function 
of  time  only  when  integrated.  A  choice  for  g(z)  which  meets  the  assumptions  is 


*(*-«.)/*  a 
s(z)  =  ^ 


(80) 


Consider  electrons  initially  located  at  height  Zq  (we  will  suppress  the  hori¬ 
zontal  dependence  of  the  electron  position;  vertical  motion  is  mainly  responsible  for 
density  changes).  The  equation  of  motion  for  the  height  of  these  electrons  at  time  t, 
z  -  z(zo,t),  in 
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z(zo,0)  =  Zq. 


(81) 


dz 

~dt  ~  i  +  w°> 

The  solution  of  this  equation  is  given  implicitly  by 


where 


z  —  ze 

2H  ’ 


and 

D  =  j  u>o(zoit>)  dt'. 


(82) 


(83) 


(84) 


The  solution  of  equation  (82),  call  it  z(zo,t ),  can  be  determined  most  readily  by  simple 
numerical  means.  The  inverse  solution,  Zo  =  Zo(z,t)  can  be  determined  in  the  same 
way  from  equation  (82),  where,  for  consistency,  z  should  replace  zo  in  equation  (84). 
The  solutions  given  are  not  exact,  but  approach  the  exact  solution  as  the  dependence 
of  w0  on  z  decreases. 

The  density  of  electrons,  initially  at  Zo ,  with  density  n0(z0) ,  are  determined 
by  dividing  by  the  jacobian  of  the  transformation,  zq  — »  z. 


n«(z,t)  = 


no(zo) 

J(zo,t)' 


(85) 


where 


J{zo,t)  = 


dz 


dzo  14 t 


l  +  <^ 

V 


(86) 


Alternatively,  the  density  can  be  written,  using  the  jacobian  of  the  inverse  transfor¬ 
mation, 
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»«(*,*)  =  n0(zo)j  (z,t), 


(87) 


where 


dz  "  l  +  ««i 


(88) 


It  can  be  readily  seen  from  these  equations  that  the  jacobians,  and  hence  the  densities, 
remain  positive  for  all  — oo  <  D  <  oo,  as  required  by  the  first  assumption.  The 
construction  also  makes  it  clear  that  for  small  velocities,  and,  hence,  small  values  of 
D,  the  solution  approaches  the  linear  solution  given  by  equation  (72),  as  demanded 
by  the  second  assumption. 

An  extension  for  N  localised  disturbances  can  be  obtained  simply  by  sum¬ 
ming  the  quantity  A  due  to  each  localized  disturbance  individually,  which  will  now 
be  called  D, 


N 

ft. 


»=1 


t— •  1,...,  N, 


(89) 


where 


A  =  f  tv0,i{ztt')dt'.  (90) 

Jo 

and  wo,i  is  vertical  velocity  of  electrons  produced  by  the  tth  localized  disturbance 
(divided  by  g).  One  then  solves  equation  (82)  to  find  the  initial  position  of  electrons 
and  equation  (87)  to  compute  the  density.  As  for  N  =  1,  all  three  of  the  required 
assumptions  are  met  by  this  solution. 
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SECTION  6 


OTHER  MECHANISMS  FOR  CHANGE  IN  ELECTRON  DENSITY 


Thus  far  we  have  been  concerned  with  changes  in  ionospheric  electron  den¬ 
sity  due  only  to  volume  changes  produced  by  motion  of  the  electrons  driven  by  the 
passage  of  an  acoustic-gravity  wave.  Electrons  can  be  freed  into  the  ionosphere  by  pho¬ 
toionization  of  atmospheric  constituents,  mainly  nitrogen  and  oxygen  (either  molecu¬ 
lar  or  atomic),  by  the  sun  (and  other,  less  important,  forms  of  radiation).  The  rate  of 
electron  production  at  a  given  location  does  not  depend  on  the  amount  of  electrons 
present,  but  rather  on  the  intensity  of  the  radiation  and  the  amount  of  atmospheric 
gases  present.  Electrons  can  be  absorbed  (or  freed)  via  photochemical  and  other  re¬ 
actions  with  atmospheric  constituents,  such  as  recombination  and  detachment.  The 
rate  of  electron  loss  due  to  such  reactions  depends  on  the  density  of  electrons  present. 
It  will  not  be  necessary  for  what  follows  to  describe  the  individual  processes  involved 
(see  reference  7  for  further  details).  In  addition  to  production  and  loss  mechanisms, 
electron  density  can  also  change  through  transport  processes,  such  as  diffusion.  Our 
object  in  this  section  is  to  amend  equation  (87)  of  the  previous  section  in  a  simple 
way  which  accounts  for  electron  density  changes  due  to  these  processes  as  modified 
by  the  passage  of  a  gravity  wave. 

The  basic  strategy  will  be  as  follows:  we  will  first  show  that  loss  mechanisms 
conspire  with  diffusion  in  such  a  way  that  the  ionosphere,  during  nighttime,  decays  in 
a  shape  preserving  manner,  which  leads  to  the  first  assumption  that  photochemical 
loss  at  all  altitudes  of  interest  depends  on  the  chemical  loss  rate  at  the  electron  density 
peak.  By  assuming  a  single  decay  rate  for  the  entire  ionospheric  profile  still  applies 
during  the  passage  of  a  gravity  wave,  we  will  then  derive  a  formula  which  gives  the 
change  in  this  rate  as  a  function  of  the  parameters  describing  the  gravity  wave.  This 
formula  will  then  be  added  to  equation  (87)  in  a  reasonable  manner.  At  this  point, 
the  reader  should  be  reminded  that  not  all  results  arrived  ?.t  in  this  section  will  be 
rigorously  derived  from  first  principles,  and  the  results  should  be  viewed  in  the  same 
spirit  as  the  previous  section. 

According  to  the  remarks  made  above,  the  continuity  equation  for  electrons 
can  be  written 


dne 

dt 


+  V  •  (n*ve)  =  q  —  I, 


(91) 
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where  the  electron  density,  ne,  and  velocity,  v4,  and  the  electron  production  rate,  g, 
are  all  functions  of  position  and  time.  The  electron  loss  rate,  l(r,t)  —  L(ne(r,f),r,t), 
is  also  a  function  of  position  and  time  through  its  dependence  of  the  function  L  on 
ne.  The  electron  velocity  can  be  thought  of  as  the  sum  of  two  velocities,  one  which 
occurs  in  ambient  conditions  due  to  transport  and  the  other  due  to  gravity  wave 
induced  electron  motion.  The  electron  production  rate  is  usually  assumed  to  be  zero 
during  the  nighttime  and  a  function  of  only  altitude  during  the  day  (a  consequence 
of  a  stratified  atmosphere),  with  a  mild  time  dependence  due  to  the  travel  of  the  sun 
across  the  daytime  sky.  A  useful  form  for  the  electron  production  rate  as  a  function 
of  altitude  is  the  classical  Chapman  formula,  which  can  be  written 


q  =  9m  exp 


1 


Z~  Zc 

H 


sec  x  e 


(92) 


where  x  is  the  solar  zenith  angle,  and  qm  is  the  production  rate  at  the  altitude  of 
the  electron  density  peak,  ze,  for  vertical  solar  radiation,  x  —  0.  This  form  for  q  is 
derived  under  the  assumption  the  that  atmosphere  is  made  up  of  a  single  absorbing 
gas  whose  density  decreases  exponentially  with  altitude  with  constant  scale  height, 
H ,  (see  reference  7). 

For  electrons  at  F- region  heights  and  above,  to  quite  a  good  approximation, 
the  electron  density  loss  rate  is  proportional  to  the  electron  density,  with  the  pro¬ 
portionality  constant  a  function  of  altitude  only,  dependent  on  the  concentration  of 
neutral  atomic  constituents.  As  shown  in  reference  7,  the  loss  rate  can  be  written 


l  =  fint,  P  =  0o  «-*-*■>/*,  (93) 

where  0o  is  loss  rate  at  altitude  zc.  The  constant,  7,  depends  on  the  relative  scale 
heights  of  the  ionizable  gas  and  linear  loss  coefficient  and  is  equal  to  about  1.75  for 
altitudes  above  the  E  region.  If  there  were  no  transport  processes  it  can  be  seen  that 
electrons  would  be  lost  at  low  altitudes  at  exponentially  increasing  rates,  which  is  not 
the  case. 


0.1  DIFFUSION. 


We  will  now  show  how  electron  diffusion  can  be  derived  from  the  equa¬ 
tions  of  motion  in  the  simplest  possible  case  and  how,  under  the  assumption  that 
the  electron  loss  rate  is  given  by  equation  (93),  this  diffusion  controls  the  decay  of 
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the  ionospheric  profile  in  such  a  way  that  its  shape  of  the  profile  is  preserved.  We 
assume  the  ionosphere  can  be  described  by  a  two  component  charged  fluid  made  of 
ions  and  electrons  undergoing  collisions  with  themselves  and  a  neutral  background. 
The  equations  of  motion  are  then 


Dv  -  1 

mt— 1  =  mjg - V(mkTi)  +  e(E  +  v,  x  B) 

JJt  n, 

-  v„)  -  mej/„(v,  -  ve),  (94) 

-  —  V{ntkTt)  -  e(E  +  v,xB) 

Dt  ne 

-mtutn{\e  -  vn)  -  mt vei(ve  -  v*),  (95) 

where  for  electrons  me  is  the  mass,  ne  the  density,  ve  the  velocity,  and  Te  the  temper¬ 
ature,  and  similarly  for  ions,  uin  is  the  collision  rate  between  ions  and  neutrals,  uen 
the  collision  rate  between  electrons  and  neutrals,  and  i/M  the  collision  rate  between 
electrons  and  ions,  E  is  the  electric  field,  B  the  magnetic  field,  g  the  acceleration  due 
to  gravity,  k  is  Boltzmann’s  constant,  and  vn  is  the  neutral  velocity. 

The  simplest  case  is  time-independent  vertical  diffusion  in  a  vertical  mag¬ 
netic  field.  In  this  case  the  equations  of  motion  reduce  to  (all  quantities  depend  only 
on  the  vertical  coordinate,  z ) 


djnjkTj) 

dz 


—  +  nitE  -  n,m,t/jn(to,  -  wn), 


(96) 


d(ntkTe) 

dz 


-nemeg  -  nteE  -  nemeutn(we  -  wn). 


(97) 


Under  the  further  simplifying  assumptions  that  »  me,  n,  =  ne  ~  n,  w,  =  we  =  Wd 
(plasma  drift  velocity),  wn  =  0,  and  rrui/in  >  mei/en,  addition  of  the  two  equations 
gives  for  the  plasma  drift  velocity 


(98) 


With  the  further  definitions 
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(99) 


JV-jpi  +  T.) 

the  plasma  drift  velocity  becomes 


2kTp 

niig  ’ 


wd 


+ 


i  dr, 
r„  dz 


+ 


where  the  plasma  diffusion  coefficient,  D,  is 


(100) 


D  s  (101) 

rriiVin 

Simplifying  further,  we  assume  the  ion  mass  and  neutral  mass  are  equal  and  are 
distributed  with  scale  height,  J5T,  and  that  Ti  =  Tt  =  T,  which  gives 


D  =  £>0e(*'")/H,  (102) 

where  X>0  is  the  diffusion  coefficient  at  zc.  The  plasma  drift  velocity  is  finally  given  by 

+  (103) 

Notice  that  the  diffusion  coefficient  increases  exponentially  with  altitude,  while  the 
electron  loss  rate,  /?,  given  in  equation  (93),  decreases  exponentially  with  altitude. 
It  is  expected  that  the  effects  will  counteract  one  another  is  such  a  way  that  a  peak 
electron  density  will  occur  near  where  /?  and  D/H 2  are  equal.  This  analysis  can  be 
generalized  by  removing  some  of  the  simplifying  assumptions  used  in  this  derivation 
(such  as  including  an  inclined  magnetic  field  and  a  nonzero  neutral  wind). 

We  will  now  try  to  find  solutions  of  the  continuity  equation  (91),  assuming 
a  loss  term  given  by  equation  (93)  and  an  electron  velocity  given  by  equation  (103). 
Assuming  vertical  motion  only  and  taking  q  =  0  for  simplicity,  equation  (91)  becomes 

dnt  „  d  .  .  ,  . 

—  =  -pnt  -  —{newd)t  (104) 
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where  ne  is  now  considered  a  function  of  z  and  t.  Assuming  a  solution  of  the  form 


nt{z,t )  =  f i(z)e~xt, 


(105) 


and  using  equation  (103)  gives 

-sK*  +  «)]+/*' (106) 

This  equation  can  be  made  into  an  eigenvalue  problem  for  A  by  supplementing  it  with 
appropriate  boundary  conditions.  We  assume  that  the  electron  density  is  zero  at  the 
ground  and  falls  off  rapidly  enough  as  z  — »  oo  to  insure  a  denumerable  number  of 
eigenvalues.  Making  the  substitution 


puts  equation  (106)  in  standard  Sturm-Liouville  form 


d_ 

dz 


+  qn=  'ffi, 


(107) 


(108) 


with  boundary  conditions,  say, 


A*(0)  =  0, 


lim^=0, 
*—00  dz 


(109) 


and  where 


p(z)  =  D{z)c-{x~“)/iH,  q{z)  =  (}[z)e-lM-M'ViH,  r{z)  =  (110) 

It  is  well  known  that  nontrivial  solutions  to  equation  (108)  exist  for  a  denumerable 
number  of  real  eigenvalues,  A„.  Multiplying  equation  (108)  by  /x,  then  integrating 
from  0  to  oo,  and  using  the  boundary  conditions  (109)  shows  that 
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>0, 


(111) 


A  = 


+ 


rp2  dz 


dz 


which  follows  from  equations  (110),  (93),  and  (102),  where  it  is  obvious  that  p,  g,  and 
r  are  positive.  Therefore,  all  eigenvalues  are  greater  than  zero.  The  general  solution  to 
equation  (104)  can  be  written  as  a  sum  of  terms  of  the  form  (105)  with  all  terms,  and 
therefore,  ne,  decaying  in  time.  Furthermore,  eventually  the  term  with  the  smallest 
eigenvalue,  A0,  will  dominate  the  sum  and  the  solution  will  rapidly  approach  the  shape 
preserving  form 


nt(z,t)  =  n0(z)e  Xot,  (112) 

It  can  further  be  shown  (see  reference  7  and  references  contained  therein)  that  the 
peak  of  the  function  n0  occurs  near  where  (3  =  D/H2,  which  we  have  denoted  /?0, 
and  that  the  smallest  decay  rate  is  approximately  equal  to  the  electron  loss  rate  at 
the  altitude  of  the  peak,  i.e.  A0  =  /?<>•  It  can  also  be  shown  that  during  the  daytime, 
q  ^  0,  an  equilibrium  solution  can  be  found  such  that  the  electron  peak  occurs  at 
approximately  the  same  altitude  and  that  the  peak  electron  density,  nm  =  n0{ze),  is 
approximately  equal  to  flo/qo ,  where  q0  is  the  value  of  the  solar  production  rate  at  the 
altitude  of  the  peak,  which  by  equation  (92)  is  q0  =  qmexp(l  —  sec  x)  for  a  Chapman 
layer. 


0.2  CHANGE  IN  LOSS  RATE  DUE  TO  GRAVITY  WAVES. 


We  have  seen  above  that  in  the  absence  of  gravity  waves  an  ionospheric 
profile  approximately  keeps  its  shape  and  decays  (in  the  absence  of  solar  radiation) 
with  a  rate  approximately  equal  to  the  electron  loss  rate  at  the  altitude  of  the  electron 
density  peak.  In  the  presence  of  solar  radiation  the  peak  approaches  a  constant  value 
approximately  equal  to  Po/go,  with  transient  terms  decaying  with  a  rate  approximately 
equal  to  Pq.  We  now  assume  that  during  and  following  the  passage  of  one  or  more 
gravity  waves,  the  electron  density  profile  is  changed  due  to  hydrodynamic  motions 
described  earlier,  but  in  addition,  the  electron  density  returns  to  its  ambient  state  in 
such  a  way  that  the  induced  transient  terms  decay  with  a  single  rate  for  the  entire 
profile,  which  is  given  by  the  ambient  rate  plus  a  correction  term  which  depends  on 
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parameters  describing  the  gravity  wave.  We  will  now  compute  an  approximation  for 
this  correction  term. 

Let  the  density  of  electrons  at  the  peak  be  denoted  nm,  and  let  the  motion 
of  the  electrons  in  the  peak  induced  by  the  gravity  waves  be  z(t),  with  z(0)  =  zc.  The 
continuity  equation  gives  for  these  electrons 

=  -£„„  +  « (log  J),  (113) 

where  all  quantities  in  this  equation  are  functions  of  time  only,  that  is, 


0{t)  =  , 


(114) 


Q{t)  =  *(*(«)),  (115) 

with  q(z)  given  by  equation  (92),  and 

dz 

m  =  (lie) 

It  will  be  assumed  that  J  is  given  by  equation  (86),  i.e.  no  other  electron  motion  of 
the  peak  is  induced  by  the  passage  of  the  gravity  wave  other  than  the  direct  coupling 
to  the  neutral  motion  given  by  equation  (77). 

The  solution  of  equation  (113)  is  readily  determined, 


_  «-*(«) 


/V<‘V(0Q(0  df  +  nmfi 
.J  0 


(117) 


where 


K(t)  =  [' fH‘‘)  dt\  (118) 

Jo 

and  rim0  =  nm(0),  the  electron  density  at  the  peak  immediately  following  the  passage 
of  the  gravity  waves,  which  is  taken  as  t  =  0. 
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For  no  gravity  wave  motion,  z(t)  =  zt ,  J(t)  =  1,  the  solution  (117)  becomes 

+S-  (119) 

It  is  clear  from  this  equation  that  if  the  initial  electron  density,  nm  0,  is  not  equal  to 
the  equilibrium  value,  ?o/A>>  then  the  solution  will  approach  the  equilibrium  value 
with  the  transient  terms  decaying  with  a  rate  equal  to  fa. 

Evaluation  of  the  integrals  appearing  in  equations  (117)  and  (118)  for  ar¬ 
bitrary  motion  can  be  performed  in  a  straightforward  manner  by  numerical  means, 
however,  this  would  not  suit  the  purposes  of  a  model.  We  therefore  seek  an  approxi¬ 
mate  solution,  which  can  readily  be  compared  with  exact  evaluation  of  the  integrals. 
The  sought  after  approximate  solution  will  be  determined  in  stages,  by  considering 
increasingly  general  cases. 

First  consider  an  undamped  gravity  wave  motion  in  tie  absence  of  solar 
radiation,  given  by  z(t)  =  zc  +  asinut,  J(t)  =  1,  Q(t)  =  0.  The  solution  (117) 
becomes 


nm  =  nmfie  K(i),  K(t)  =  0QI(t),  (120) 


where 


I(t)  =  J*  dt\ 

and  b  =  7a.  It  can  easily  be  shown  that 


(121) 


J(«)  =/„(»)(! -r)  +  oscillatory  terms,  0  <  ur  <  2ir,  (122) 

where  la  is  the  zeroth  order  modified  Bessel  function  of  the  first  kind.  For  small  6, 
the  power  series  expansion  of  Jo  gives. 


(123) 
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so  that  the  solution  (120)  becomes 


nm  =  nmio  e’*t  0  =  fa  (l  +  ,  (124) 

Notice  that  for  a  single  undamped  gravity  wave  with  amplitude  not  too  large,  the 
peak  electron  density  decays  with  a  rate  equal  to  the  ambient  rate  plus  a  second  order 
correction  that  depends  only  on  the  amplitude  of  the  wave  and  not  on  the  frequency. 
There  is  no  first  order  correction,  which  shows  that  the  amplitude  must  be  somewhat 
large  to  see  a  change  in  decay  rate. 

An  alternative  way  of  calculating  the  correction  to  the  decay  rate  is  to 
expand  the  integrand  in  equation  (121) 

I(t)  =  /  Y*  (— fcsinwf')"  dt '.  (I25) 

Jo  *=o 

Retaining  terms  to  second  order  gives 

I(t)  =  ^1  +  ~  j  cos  wt  ^1  —  ^  sinwf^  “  ~  - >  (126) 

which  gives  the  same  second  order  change  in  the  decay  rate  as  in  equation  (123). 

We  now  generalize  this  result  by  considering  a  damped  gravity  wave,  such 
that  z(t)  —  ze  +  ae~xt  sin  wt,  J(t)  =  1,  Q(t)  —  0.  The  solution  is  again  given  by 
equation  (120)  with 


J(i)  =  J  exp  j— be  xt'  sinwt'j  dt', 
Expanding  in  the  same  way  as  equation  (125)  gives 


(127) 


(128) 
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Notice  that  this  gives  the  same  result  as  in  equation  (126)  in  the  limit  as  X  — ►  0.  This 
solution  gives  a  decay  rate  which  initially  has  a  second  order  increase,  but  eventually 
returns  to  its  ambient  value,  as  would  be  expected  as  the  gravity  wave  damps  out. 

Generalizing  further  let  the  electron  motion  be  given  by  the  sum  of  N 
damped  gravity  waves  with  differing  frequencies  and  phases, 

AT 

z{t)  =  Ze  +  Yi,  8in(w,-t  +  t  =  1, . . . ,  N.  (129) 

«=i 


The  solution  is 


rt  r  n 

I(t)  =  j  exp  |  —  Y  sin(a>1t*  +  <f>j)  dt‘,  (130) 

With  a  bit  more  algebra  one  can  show  that  to  second  order,  for  u,  ^  u>j,  i  ^  j, 


m 


i  n 


b} 


w; 


8  ,=i  K  + 


(131) 


Most  interesting  about  this  result  is  that  if  the  frequencies  of  oscillation  are  different, 
then  the  second  order  correction  is  given  by  the  sum  of  terms  which  apply  for  one 
oscillation  (as  in  equation  (128)),  there  is  no  mixing  of  terms  (no  terms  in  6,6,  ,  x  ±  j). 

If  the  gravity  wave  is  damped  using  a  gaussian,  as  was  chosen  in  equa¬ 
tion  (74), 


N 


z{t)  -  ze  +  ^2<ne~^1  sin(u/,f  +  <£,), 

»=i 


*  =  1,... ,  N, 


(132) 


then,  assuming  A,  <  w,, 


/(«)=<+ I 


(133) 


This  is  the  most  general  form  we  will  use  for  the  change  in  the  decay  rate  due  to  the 
presence  of  many  damped  gravity  waves. 
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Performing  the  more  general  calculations  when  Q(t)  ^  0  and  J(t)  ^  1 
becomes  prohibitively  more  difficult.  It  tur«.s  out  that  to  the  order  in  which  we  are 
working  we  can  approximate  the  solution  (117),  by 


nm  = 


(134) 


where 


fi{t)  =  (30I(t)/t.  (135) 

A  comparison  of  the  approximate  formula  (134)  with  the  exact  result  (117), 
is  shown  in  figure  7.  The  electron  density  has  been  converted  to  critical  frequency 
using  equation  (1)  for  a  case  similar  to  that  shown  in  figure  1.  As  can  be  seen,  the 
approximate  result  (dashed  line),  agrees  well  with  the  exact  integral  (solid  line). 
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Figure  7 
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,  A  comparison  of  critical  frequency  vs.  time  using  the  approximate 
formula  (134)  with  exact  integral  (117). 
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SECTION  7 


THE  MODEL 


We  will  now  assemble  all  of  the  ideas  derived  in  the  previous  sections  into  a 
working  model.  It  will  be  necessary  for  both  consistency  and  simplicity  to  make  some 
changes  in  the  formulas  presented  earlier,  as  well  as  to  include  features  not  previously 
discussed.  All  the  derivations  were  carried  out  in  cartesian  coordinates,  assuming  a 
flat  earth,  with  an  atmosphere  of  constant  temperature  and,  therefore,  constant  scale 
height  and  speed  of  sound.  The  model  must  compute  electron  densities  above  a  round 
earth,  so  the  necessary  geometrical  changes  will  be  made.  The  actual  atmosphere  has 
a  scale  height  which  increases  with  altitude.  For  simplicity,  we  will  choose  an  average 
scale  height  and  speed  of  sound,  which  give  results  that  best  approximate  the  data. 
In  section  4.4  it  was  shown  that  the  amplitude  of  gravity  wave-induced  ionospheric 
disturbance  falls  off  as  the  inverse  distance  from  the  source  under  the  assumptions 
made  in  that  section.  In  a  more  realistic  case,  there  will  be  other  mechanisms  which 
will  cause  damping  of  the  wave.  We  will  add  an  exponential  damping  term  to  account 
for  this  expected  behavior,  with  the  length  scale  chosen  to  best  agree  with  the  data. 

The  ambient  electron  density  profile  for  the  ionosphere  in  general  depends 
on  location  and  time  of  day  (as  well  as  other  factors  such  as  season,  solar  activity,  etc.). 
For  modeling  purposes  we  will  use  two  vertically-stratified  electron  density  profiles  to 
represent  typical  day  and  nighttime  cases.  The  altitude  of  the  -lectron  density  peak, 
ze,  will  be  taken  to  be  300  km  and  the  ambient  electron  loss  rate  at  this  altitude, 
/?o,  will  be  chosen  to  reproduce  typical  nighttime  behavior.  The  daytime  electron 
production  rate  qo  will  be  chosen  such  that  no  =  qo/Po,  where  n0  is  daytime  electron 
density  profile. 


7.1  CONSTANTS  AND  AUXILIARY  FUNCTIONS. 


The  following  list  of  constants  have  been  chosen  to  best  fit  the  available  data 
(they  may  be  changed  by  the  user): 


c  —  .7  km/sec,  ze  =  300  km,  H  =  40  km,  (136) 
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x0  —  1  x  104  km, 


A0  =  50  km/MT,  m0  =  .45  MT"1/2,  (137) 


/J0=lx  10-4  sec-1,  7  =  1.75.  (138) 

It  is  assumed  that  the  model  is  provided  with  a  number  of  external  functions 
which  can  be  called  when  necessary.  They  are: 

Electron  density  profile:  no (z)  (one  for  day  and  one  for  night), 

Brunt- Vaisala  frequency:  wj(a). 

Acoustic  cutoff  frequency:  ua(z), 

A  function:  zq(z,D)  and  its  inverse  z(zq,  D),  which  are  solutions  of 


where 


z'  —  t  * 


~  e 


-*o 


+  D', 


Z-Zc  ,  _  Zp-  Zc 

H  ’  JT  ’ 


(139) 


(140) 


7.2  INPUTS  AND  OUTPUTS. 


The  earth  is  taken  to  be  a  sphere  of  radius  Re,  with  a  spherical  coordinate 
system,  {r,0,  <f>),  centered  at  the  center  of  the  sphere.  The  earth’s  magnetic  field  is 
taken  to  be  a  magnetic  dipole  with  axis  along  the  polar  axis  of  the  spherical  coordinate 
system. 


The  inputs  to  the  model  are  the  locations,  yields,  and  times  of  burst  of 
N  explosions,  and  the  location  and  time  at  which  the  electron  density  is  desired 
(observation  point).  The  following  list  describes  each  input  parameter: 

0,,  i  =  1, . . . ,  N:  colatitude  of  the  »th  burst, 

4>i,  t  =  1, . . . ,  N:  longitude  of  the  *th  burst, 
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zhji,  i=  1  altitude  of  the  *th  burst  in  km, 

Yi,  i  =  1, . . . ,  N:  yield  of  the  *th  burst  in  MT, 
to,,*,  i  =  1  time  of  the  »th  burst  in  seconds, 

6:  colatitude  of  the  observation  point, 

<f>:  longitude  of  the  observation  point, 
z:  altitude  of  the  observation  point, 
t:  time  of  observation. 

The  outputs  of  the  model  are  the  electron  density  as  a  function  of  position 

and  time, 


MM,  <t> ,  0, 

each  partial  derivative  of  ne  with  respect  to  r,  0,  <j>,  and  t, 

dne  dne  dnt  dne 

H7'  1$'  Ht'  Ht' 

as  well  as  all  second  order  partial  derivatives  of  nt  with  respect  to  r,  6,  and  <f>, 

d2nt  d2ne  d3ne  d?ne  d2ne  d2ne 

Hr2'  HP'  Htf'  IfrdO'  dr<f> '  dOdf 

We  will  only  present  the  formulas  for  the  electron  density,  ne,  the  formulas  for  the 
first  and  second  partial  derivatives  can  be  obtained  from  these  formulas  by  a  straight¬ 
forward,  though  tedious,  application  of  the  chain  rule  of  the  differential  calculus  of 
many  variables. 


7.3  MODEL  ALGORITHM. 

The  electron  density,  MM,<M)»  is  computed  as  follows: 

rh,i  ~  Re  t  ,  *  !,•••>  (141) 
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Z{  —  Z  Zb,iy 

(142) 

cos  $Zii  —  cos  6  cos  Oi  4-  sin 9  sin  0,-  cos(<f>  —  <f>i) , 

(143) 

= 

(144) 

iZ?  =  xf  +  2-, 

(145) 

.  x  2  cos  6 

Sm  (1  +  3  cos*  tf)1/1’ 

(146) 

cos  0,  —  cos  0  cos  0S>< 
t,to0i  sin  9  sin  0*it 

(147) 

P,  =  -  (cos  a,  sin  2A-^-  —  2  sin*  , 

2  \  Ri  R%J 

(148) 

Ei  =  c~*i/x\ 

(149) 

,  Ri“4z)  . . 

•l,i  —  /  \  '  l0,»j 

C  U>4  (z) 

(150) 

.  «<  /  lwo(z)  Pf  ^ 

“*w *  V  +  H’W-‘(| -<o,)V’ 

(151) 

Qi  =  sin  (wi(t  -  ti,i)) , 

(152) 

Ai  -  A0Yi,  m,  =  moy/Y» 

(153) 
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27T  mi 

Oi  = 

2. 

(154) 

Ri 

ki  =  Wi(z)(f  -  «!,.)— 2, 

Xi 

(155) 

Si- Ai  ex  p  (~k]ja})  , 

(156) 

A  =  S,QiEiPi, 

(157) 

N 

c  =  E  A, 

«=1 

(158) 

di  =  A<APi, 

(159) 

^  Zc  24  ,- 

*  2jt  m,-  x,-  ’ 

(160) 

<16» 

&  =  +%/&*, 

(162) 

/?  =  &  |l  +  [z(ze,  -6)  -  *<])  ] 

,  (163) 

/  A>n0(2),  day 

90  ”  \  0,  night  » 

(164) 

Zq  =  Zq{z,D), 

(165) 
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(166) 


T  1  +  e~< 

**  —  „  J 

1  +  c 

i  =  t  —  min^ii),  (167) 

ne  =  ^no(zo)  “  j)  «"*  +  J  j  J-  (168) 

7.4  RESULTS. 

Figures  8  and  9  are  electron  density  (plasma  frequency)  contour  plots  gen¬ 
erated  using  the  AGW  model.  For  each  plot,  the  electron  density  following  the 
detonation  of  four  simultaneous  high-yield,  near  surface  bursts  at  two  hours  is  plotted 
for  the  surrounding  region.  Figure  8  is  the  daytime  case  and  figure  9  is  the  nighttime 
case. 

The  lower  portions  of  each  figure  are  electron  density  contours  in  a  horizontal 
plane  300  km  above  the  earth’s  surface.  The  distance  between  the  tick  marks  is 
556  km.  The  upper  portions  of  each  figure  are  electron  density  contours  in  a  vertical 
plane  which  intersects  the  lower  portion  of  the  figure  along  the  line  shown  in  the  center 
of  the  figure.  It  is  clear  from  the  figures  that  the  effect  of  the  gravity  waves  is  to  distort 
the  normally  vertically  stratified  ionosphere  into  regions  of  varying  density  on  scales 
of  about  10  km  to  100  km.  Notice  that  the  effect  is  more  prominent  during  nighttime. 
During  the  daytime,  the  sun  acts  to  more  swiftly  return  the  ionosphere  to  its  ambient 
state.  Further  details  concerning  these  plots,  as  well  as  the  implementation  and  use  of 
the  AGW  model  in  the  RAYTRACE  code  and  its  subsequent  application  for  an  HF 
signal  specification  can  be  found  in  the  companion  report,  reference  6. as 
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SECTION  9 


GLOSSARY  OF  SYMBOLS 


a  2  r mfz  and  AGW  amplitude 

b  7  a 

B  earth’s  magnetic  field 

c  speed  of  light,  or  sound 

cv  specific  heat  (constant  volume) 

e  (w*/w0)  c 

D  plasma  diffusion  constant 

D0  plasma  diffusion  constant  (reference) 

f  wave  frequency 

foF2  critical  frequency  of  F2  layer 

ip  plasma  frequency 

f(r)  initial  velocity  distribution 

/( k)  fourier  transform  of  f(r) 

g  acceleration  of  gravity 

g(z)  height  function 

h  true  height 

H  scale  height 

h'  virtue  height 

1(A)  stationary  phase  integral 

j  inverse  jacobian 

J  jacobian 

k  wave  number 

K  constant 

k,  ith  component  of  k 

k*  horizontal  wave  number 

k0  stationary  phase  solution 

l  electron  loss  term 

m  number  of  cycles 

me  electron  mass 

m,  ion  mass 

m0  number  of  cycles  at  yield  one 

n  index  of  refraction 

N  number  of  disturbances 

n’  perturbation  in  electron  density 


ne(z)  electron  density  at  height  z 

n e.mox  ionospheric  peak  electron  density 

n,-  ion  number  density 

nm  electron  density  at  peak 

nm,^  electron  density  at  peak  just  after  AGW 

n^  ambient  electron  density 

p  pressure 

pW  ith  pressure  perturbation 
Po  pressure  at  height  zero 

q  electron  production  rate 

q^  electron  production  rate  at  peak 

r  radius 

R*  earth’s  radius 

Th  cylindrical  radius 

s  entropy 

t  time 

Te  electron  temperature 

T,-  ion  temperature 

Tp  plasma  temperature 

T*  disturbance  onset  time 

u  horizontal  speed 

U  constant 

u^  ith  horizontal  speed  perturbation 
u,  3u fdx  (etc...) 

U„  U  at  Y  1 

v  velocity 

vt  electron  velocity 

vt  group  velocity 

Vi  ion  velocity 

v(’)  ith  vertical  speed  perturbation 
w  vertical  speed 

Wd  plasma  drift  velocity 

x  horizontal  position 

Y  yield 

z  height 
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Zc 

P 

P* 

1 

V{z) 

9 

£ 

An 

A0 

A 

I'in 

P 


X 

T 

We 
Wj 
W  a 
W 
W  „ 
Wc 


scaled  height 
reference  height 
t/ne 

loss  rate  reference 
ratio  of  specific  heats 
functional  form  of  n«(z) 
angle  w.r.t.  vertical 
expansion  parameter 
eigenvalues  of  A 
dominant  eigenvalue 
parameter 

normalized  functional  form 
ion  neutral  collision  frequency 
density 

density  at  height  zero 
ith  density  perturbation 
zenith  angle 
ct/r 

vertical  velocity  of  electrons 
Brunt- Vaisala  frequency 
acoustic  cutoff  frequency 
angular  frequency 
stationary  phase  frequency 
W&  cos  9 
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